library(dplyr)
library(magrittr)
library(tidyr)
library(ggplot2)
library(ggmap)
library(maps)
source("vcode_function.R")

load("data/demography_data_smoothed.Rdata")

eur <- dat %>%
  select(country_id, eur_pct_est_smooth, year)

world <- map_data("world")
world %<>%
  filter(!region %in% c("French Southern and Antarctic Lands",
                        "Siachen Glacier",
                        "Antarctica"))
world <- vcode(world, "region", custom_matches = c("Virgin Islands" = "U.S. Virgin Islands"))

epes_1600 <- eur %>%
  filter(!is.na(eur_pct_est_smooth) & year == 1600) %>%
  group_by(country_id) %>%
  do(slice(., 1)) %>%
  ungroup %>%
  mutate(country_id = as.numeric(country_id)) %>%
  select(-year)
w0 <- left_join(world, epes_1600)
world_map_1600 <- ggplot(data = w0,
                         mapping = aes(x = long,
                                       y = lat,
                                       group = group)) + 
  geom_polygon(aes(fill = eur_pct_est_smooth), color = "black") +
  scale_fill_gradient(name = "European Ancestry", na.value = "white", low = "white", high = "black") + 
  coord_fixed(1.3) +
  labs(title = "1600", x = NULL, y = NULL) +
  theme(legend.position = "bottom",
        panel.background = element_rect(color = "black", fill = "white"),
        axis.text = element_blank(),
        axis.ticks = element_blank())

epes_1700 <- eur %>%
  filter(!is.na(eur_pct_est_smooth) & year == 1700) %>%
  group_by(country_id) %>%
  do(slice(., 1)) %>%
  ungroup %>%
  mutate(country_id = as.numeric(country_id)) %>%
  select(-year)
w1 <- left_join(world, epes_1700)
world_map_1700 <- ggplot(data = w1,
                         mapping = aes(x = long,
                                       y = lat,
                                       group = group)) + 
  geom_polygon(aes(fill = eur_pct_est_smooth), color = "black") +
  scale_fill_gradient(name = "European Ancestry", na.value = "white", low = "white", high = "black") + 
  coord_fixed(1.3) +
  labs(title = "1700", x = NULL, y = NULL) +
  theme(legend.position = "bottom",
        panel.background = element_rect(color = "black", fill = "white"),
        axis.text = element_blank(),
        axis.ticks = element_blank())

epes_1800 <- eur %>%
  filter(!is.na(eur_pct_est_smooth) & year == 1800) %>%
  group_by(country_id) %>%
  do(slice(., 1)) %>%
  ungroup %>%
  mutate(country_id = as.numeric(country_id)) %>%
  select(-year)
w2 <- left_join(world, epes_1800)
world_map_1800 <- ggplot(data = w2,
                         mapping = aes(x = long,
                                       y = lat,
                                       group = group)) + 
  geom_polygon(aes(fill = eur_pct_est_smooth), color = "black") +
  scale_fill_gradient(name = "European Ancestry", na.value = "white", low = "white", high = "black") + 
  coord_fixed(1.3) +
  labs(title = "1800", x = NULL, y = NULL) +
  theme(legend.position = "bottom",
        panel.background = element_rect(color = "black", fill = "white"),
        axis.text = element_blank(),
        axis.ticks = element_blank())

epes_1900 <- eur %>%
  filter(!is.na(eur_pct_est_smooth) & year == 1900) %>%
  group_by(country_id) %>%
  do(slice(., 1)) %>%
  ungroup %>%
  mutate(country_id = as.numeric(country_id)) %>%
  select(-year)
w3 <- left_join(world, epes_1900)
world_map_1900 <- ggplot(data = w3,
                         mapping = aes(x = long,
                                       y = lat,
                                       group = group)) + 
  geom_polygon(aes(fill = eur_pct_est_smooth), color = "black") +
  scale_fill_gradient(name = "European Ancestry", na.value = "white", low = "white", high = "black") + 
  coord_fixed(1.3) +
  labs(title = "1900", x = NULL, y = NULL) +
  theme(legend.position = "bottom",
        panel.background = element_rect(color = "black", fill = "white"),
        axis.text = element_blank(),
        axis.ticks = element_blank())

epes_1950 <- eur %>%
  filter(!is.na(eur_pct_est_smooth) & year == 1950) %>%
  group_by(country_id) %>%
  do(slice(., 1)) %>%
  ungroup %>%
  mutate(country_id = as.numeric(country_id)) %>%
  select(-year)
w5 <- left_join(world, epes_1950)
world_map_1950 <- ggplot(data = w5,
                         mapping = aes(x = long,
                                       y = lat,
                                       group = group)) + 
  geom_polygon(aes(fill = eur_pct_est_smooth), color = "black") +
  scale_fill_gradient(name = "European Ancestry", na.value = "white", low = "white", high = "black") + 
  coord_fixed(1.3) +
  labs(title = "1950", x = NULL, y = NULL) +
  theme(legend.position = "bottom",
        panel.background = element_rect(color = "black", fill = "white"),
        axis.text = element_blank(),
        axis.ticks = element_blank())

epes_2000 <- eur %>%
  filter(!is.na(eur_pct_est_smooth) & year == 2000) %>%
  group_by(country_id) %>%
  do(slice(., 1)) %>%
  ungroup %>%
  mutate(country_id = as.numeric(country_id)) %>%
  select(-year)
w4 <- left_join(world, epes_2000)
world_map_2000 <- ggplot(data = w4,
                         mapping = aes(x = long,
                                       y = lat,
                                       group = group)) + 
  geom_polygon(aes(fill = eur_pct_est_smooth), color = "black") +
  scale_fill_gradient(name = "European Ancestry", na.value = "white", low = "white", high = "black") + 
  coord_fixed(1.3) +
  labs(title = "2000", x = NULL, y = NULL) +
  theme(legend.position = "bottom",
        panel.background = element_rect(color = "black", fill = "white"),
        axis.text = element_blank(),
        axis.ticks = element_blank())


my_grid <- ggpubr::ggarrange(world_map_1600,
                             world_map_1700,
                             world_map_1800,
                             world_map_1900, 
                             world_map_1950, 
                             world_map_2000,
                             ncol = 2,
                             nrow = 3,
                             common.legend = T,
                             legend = "bottom")

## ggsave(my_grid,
##        file = "figure_9_2.png",
##        height = 11,
##        width = 8.5)

ggsave(my_grid,
       file = "output/map_9_1.tiff",
       height = 11,
       width = 8.5,
       dpi = 300)
